Diffusion and Relaxation Dynamics in Cluster Crystals 
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For a large class of fluids exhibiting ultrasoft bounded pair potentials, particles form crystals 
consisting of clusters located in the lattice sites, with a density-independent lattice constant. Here 
we present an investigation on the dynamic features of a representative example of this class. It is 
found that particles can diffuse between lattice sites, maintaining the lattice structure, through an 
activated hopping mechanism. This feature yields finite values for the diffusivity and full relaxation 
of density correlation functions. Simulations suggest the existence of a localization transition which 
is avoided by hopping, and a dynamic decoupling between self- and collective correlations. 

PACS numbers: 61.20.Ja, 82.70.-y, 64.70.Pf 



The investigation of large-scale structural and dy- 
namic properties of complex fluids can be facili- 
tated by coarse-graining intramolecular fast degrees- 
of-freedom. By following this procedure, each macro- 
molecule is represented as a single particle interact- 
ing with any other through an effective ultrasoft pair 
potential pj]. For isotropic interactions, the latter 
just depends on the distance between centers-of-mass. 
Details of macromolecular size and strength of the 
monomer-monomer interactions enter into parame- 
ters of the effective potential. The latter is bounded 
if centers-of-mass can coincide without violating ex- 
cluded volume conditions. Some examples are poly- 
mer chains Q , dendrimers Q , or microgels 0] . 

Generalized exponential models (GEM), v{r) = 
eexp[— (r/cr)™], constitute a class of such effective 
bounded interactions. The cases m < 2 and m > 2 be- 
long, respectively, to the so-called and classes 
for which the Fourier transform, v{q), of v{r) is, re- 
spectively, positive definite or oscillating around zero. 
According to a general criterion based on a mean- 
field density functional theory systems belonging 
to the (5+ class display reentrant crystallization in the 
density-temperature plane [il, (tI. Potentials of the 
class do not yield reentrance but rather a monotonic 
freezing line beyond which the system forms cluster 
crystals fliSll^. These are novel forms for the self- 
organization of (soft) matter [21 in the sense that they 
feature a lattice constant that is density-independent 
and the cluster population is proportional to the den- 
sity 0, 12, . Recent simulations of amphiphilic den- 
drimers show that the latter indeed interact via a Q^- 
potential [l^ . 

Previous theoretical and computational investiga- 
tions of the GEM model for to > 2 have focused on 
structural and thermodynamic properties [H, [l^ • The 
dynamics of these crystals have not been studied to- 
date and their particular form of self-organization sug- 
gests that the former must be different from that of 
atomic solids with single site occupancy. The purpose 
of this Letter is to investigate the novel features aris- 
ing in the dynamics of such cluster crystals in gen- 
eral. We find that particles diffuse between lattice 
sites through an activated hopping mechanism, with- 
out breaking the lattice structure, resulting in finite 
values for the diffusivity and full relaxation of self- 



correlation functions. We establish the existence of 
a localization transition which is avoided by hopping 
and is distinct for self- and density-density correla- 
tions. Therefore, an unusual decoupling between self- 
and collective dynamics is observed, the latter exhibit- 
ing a lower transition temperature than the former. 

Particles in the simulated system interact through 
a GEM potential with to = 8, which is cut-off at a 
distance rc = 1.5a. The particle mass mp, energy 
scale e, and particle diameter a are set to one. In 
the following, energy, temperature, density, time, dis- 
tance, and wavevector will be given respectively in 
units of e, e/fce, a~^, cr(TOp/e)^/^, a, and a~^. We 
investigate dynamic features of the cluster crystals, 
in a wide range of temperature T, deep inside the 
region of stability of the fcc-phase P, [loj . The inves- 
tigated densities are p — N/L^ — 2.0, 3.0, 4.0, and 
7.0. The values of the cell size, L, used for periodic 
boundary conditions are integer multiples of the p- 
independent lattice constant a. Following d, for 
a fee lattice a = 27r\/3/9*, with the wavevector 
for which v{q) shows its first minimum. The obtained 
value is a = 1.8571325. The distance between nearest- 
neighbor lattice sites is d„„ a/%/2 = 1.313191. 

Due to the large values of p, investigation of relax- 
ation dynamics in the present system is computation- 
ally more demanding that for standard model systems 
with unbounded interactions as, e.g., Lennard Jones 
mixtures [l_3|i which are typically studied at p ^ 1.0. 
For this reason Newtonian dynamics simulations have 
been carried out instead of Brownian dynamics, which 
would yield a much slower relaxation 14 1 . 

Fig. [T] shows the T-dependence of the mean squared 
displacement for densities p = 2.0 and p = 7.0. 
The initial ballistic regime corresponds to oscillations 
around lattice sites whose period scales, according to 
[loj as p~^l^. Rescaling time as p^/^i causes a col- 
lapse of the curves at short times and confirms the 
above result that was derived on the basis of statics 
[lo| . Short-time motion is dictated, thus, by p alone 
whereas long-time diffusion by p/T, as we will shortly 
demonstrate. After the ballistic regime, a plateau 
arises characterizing the temporary trapping of the 
particles within the clusters. At long times, the par- 
ticles reach the diffusive regime ((Ar^(i)) ~ t) at all 
the investigated temperatures, reaching distances of at 
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FIG. 1: Temperature dependence of the mean squared 
displacement for density p — 2.0 (symbols) and 7.0 (lines). 
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FIG. 3: Projection of the potential energy [7 on a [001] 
lattice plane, for a configuration at p = 5.0 and T = 0.667. 
Darker colors represent lower energies (see legend) and are 
centered around lattice sites. 
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FIG. 2: The van Hove self-correlation function, dnnGs (r, t) 
computed at different times, vs. distance r/dnn, where dim 
is the distance between nearest-neighbor lattice sites, at 
p = 2.0 and T = 0.267. Arrows correspond to distances 
between nth-neighbor lattice sites, for n = 1, 2, 3, and 4. 



FIG. 4: The distribution of local potential energy, pg(U), 
vs. rescaled energy, U/p, at different densities, (a): results 
considering all the volume in the simulation box; (b): re- 
sults by computing U only at lattice sites. 



least one particle diameter from their initial position. 
The fee structure remains stable in the whole time 
window of the simulation. Therefore, particles move 
between neighboring clusters, modifying thus the clus- 
ters' initial identity but leaving unaffected their aver- 
age population Uc- 

Fig. [2] displays an illustrative example of the time 
evolution of the van Hove self-correlation function, 
Gs (r, t) . The discrete nature of the motion between 
clusters centered around distinct lattice sites is evi- 
denced by the succession of peaks arising in Gs{r,t). 
With increasing time, the height of the first peak pro- 
gressively decreases and new peaks of increasing inten- 
sity located at larger distances arise. From integration 
of Gs{r, t) for the longest represented t it is found that 
more than a 80 % of the particles are located at dis- 
tance r > O.Sdnn from their initial position, i.e., they 
have moved to different lattice sites. 

Though a detailed characterization of jumps be- 
tween lattice sites is beyond the scope of this Letter, 
it is worth commenting some interesting features of 
Gs (r, t) . Arrows in Fig. [2] indicate distances between 
?ith nearest-neighbor lattice sites, ^/nd^-a_, where n = 
1, 2, 3, and 4. The two first sharp maxima after the 
initial peak match such distances for n = 1 and 3. The 
pronounced minimum at r = 1.4dnn matches the case 
n = 2. Given two consecutive jumps connecting three 
distinct lattice sites 1, 2, 3 for which di2 = d23 = dnn, 



it is easy to see that the four allowed angles between 
both jumps, 60°, 90°, 120°, and 180° yield, respec- 
tively, dia/dnn = 1, V2, VS, and 2. The presence 
of a marked minimum for r/dnn = and a sharp 
maximum for r/dnn = suggests a preferential di- 
rectionality for the motion between neighboring sites 
with low and high probability for angles of, respec- 
tively, 90° and 120° between consecutive jumps. 

Fig. [3] shows for a typical configuration at p = 5.0 
and T — 0.667 a map of the local potential energy 
[/ at a [001] lattice plane, created by all the parti- 
cles in the system. According to Fragner [l^, for a 
sufficiently small lattice constant the superposition of 
potentials of the class centered at crystal lattice 
sites exhibits potential minima and maxima at respec- 
tively the lattice and interstitial sites. Results in Fig. 
[3] confirm this prediction. The potential barriers be- 
tween nearest-neighbor minima (separation dnn) are 
much lower than between next-nearest ones (separa- 
tion a) , suggesting a much lower probability for single 
jumps of length a. This feature is consistent with the 
minimum at \/2c?nn = ci observed in Gg {r, t) for all the 
densities (see Fig. [5]) 

Fig. m shows the distribution of local [/-values, com- 
puted over the whole simulation box [gt{U), panel(a)] 
and by only considering lattice sites [g\{U), panel(b)]. 
The distribution g\{U) is dominated by the contri- 
bution of the individual cluster centered around the 
considered site. Since the number of particles form- 
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FIG. 5: Temperature dependence of the difFusivity for dif- 
ferent densities. Lines are fits to Eq. ([T]); the values of the 
fit parameter 70 are indicated on the plot. 



ing the cluster is polydisperse and obviously discrete, 
gi{U) shows a peak structure. Consistently with the 
common value of the lattice constant and the rela- 
tion Uc oc p, approximate p-scaling behavior is ob- 
served for gt{U), of increasing quality at higher densi- 
ties. The maximum at Ut = 4/3 corresponds to typical 
values at interstitials between nearest-neighbor sites 
(see Fig. Though the most probable U for gi{U) 
shows a weak p-dependence, we take a 'common' value 
?7i = 1.7p as a rough estimation of the average local 
minimum at the lattice sites. Hence, the typical po- 
tential barrier separating nearest-neighbor lattice sites 
is AU = Ut-Ui^ 2.3p. 

Fig- m displays results for the diffusivity D, obtained 
as the long-time limit of (Ar^(t))/6i. The depen- 
dence of D on density and temperature can be un- 
derstood as follows. Each particle located on any 
given lattice site experiences a harmonic potential 
with a barrier height AU to the neighboring site, as 
mentioned above. Due to the contact with the heat 
bath, the energy i? of a particle follows the distri- 
bution f{E) = e~'^'^{f3^E'^)/2 for £; > and van- 
ishing for E < 0, with /? = (kBT)-^. The proba- 
bility P>(A[/) to have E > AU can be calculated 
as P>iAU) = e~^^^^[{/3AU)^/2 + pAU + 1]. Once 
E > AU, the particle hops to the neighboring site. 
We can now translate P>(AC/) into a 'waiting time' 
that particles typically spend on a given lattice 
site before hopping to the next one, which scales as 
t-w ^ ^/P>{AU). If we observe any given particle for 
a number of timesteps Ng ^ 1, there will only be 
a fraction A^hop ~ P>iAU)Ns of them for which the 
particle will hop, its motion being a random walk of 
step size i ^ a. Hence, measured in the natural units 
of the problem, {Ar^{t)) will scale as P>{AU)t. By 
using AU = 2.3/9 (see above) we obtain 

£» = 70 [{2.3p/Tf/2 + 2.ip/T + l] e^^-Sp/T^ 

with some numerical coefficient 70 of order unity. The- 
ory predicts, thus, that D depends solely on the ratio 
p/T. The approximations involved in deriving Eq. ([T]) 
above become more accurate as p grows, for which 
case the polydispersity in the cluster population is re- 
duced and all sites can be treated as identical har- 
monic wells. The results shown in Fig. [5] fully corrob- 
orate this treatment. The diffusivity can be very well 
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FIG. 6: Temperature dependence of Fs{q,t) [(a)] and 
F{q,t) [(b)], for p = 2.0 and q = 4.0. The insets show 
a magnification of the plateau regime. 



described by the law of Eq. ([T]) above and, indeed, for 
sufficiently high densities all data points collapse on a 
single curve, where D is expressed as a function of the 
ratio p/T. This is a dynamical generalization of the 
scaling properties previously found for statics [lo| . 

Finally, we discuss the self- and density-density 
correlators of wavevector q, respectively defined as 
F,{q,t) = iV-i(E^.exp[iq • (r,(t) - r,(0))]), and 
F{q,t) = {pgit)PgiO))/{pMP*,m, where p,{t) is de- 
fined as pq{t) = J2j exp[iq • rj{t)]. Fig. [6] shows, for 
p = 2.0 and for a fixed wavevector q — A.O not probing 
the lattice structure, the T-dependence of Fs{q, t) and 
F{q,t). After the ffist microscopic decay, both cor- 
relators exhibit a plateau, whose duration increases 
with decreasing T. As usually found in systems dis- 
playing slow relaxation, the presence of this plateau 
indicates a temporary freezing of such correlations. 
Strong density oscillations related with intracluster 
motion are observed at the crossover between the mi- 
croscopic and plateau regimes [3| • At long times self- 
correlations relax and decay to zero. This is also true 
for density-density correlations if q does not match 
any wavevector probing the lattice structure. On the 
contrary, and consistently with the observed stabil- 
ity of the fee lattice, density-density correlations for 
q belonging to the reciprocal lattice are permanently 
frozen, and there is no signature of a final decay of 
the plateau, which remains flat (not shown) up to the 
limit of the simulation window. 

An ideal (i.e., without intervening hopping events) 
localization transition to a non-ergodic ('glassy') 
phase can be defined by a jump from zero to a fi- 
nite value (non-ergodicity parameter) of the long-time 
limit of dynamic correlators [l3|. Hence, at the transi- 
tion point the non-ergodicity parameter is equal to the 
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g-dependence of the non-ergodicity parameters 
fq, lines) and density-density symbols) cor- 
Data for /| are displayed for fixed p = 2.0 and 



several temperatures. Data for fq are shown for different 
densities and correspond to state points above the local- 
ization transition. Data for fq are only shown for g-values 
not belonging to the reciprocal lattice. 



plateau height. This quantity progressively increases 
with decreasing temperature in the non-ergodic phase. 
The presence of hopping events can, in principle, re- 
store ergodicity at temperatures below the ideal tran- 
sition, leading to final relaxation of correlators [17]. 
Still, the existence of an avoided ideal transition can 
still be identified by the increase of the plateau height. 
Now we investigate features of such a transition in 
the present system. As shown in Fig. [51 Fs{q,t) and 
F{q, t) exhibit a rather different behaviour at the in- 
termediate plateau regime. The plateau height, /|, 
for Fs{q,t) increases with decreasing T [see inset in 
panel (a)]. This feature is observed for all wavevec- 
tors and suggests the existence of a localization tran- 
sition for self-motions. Fig.[7]shows the q-dependence 
of for p = 2.0 at the same temperatures of Fig. [6l 
The increasing width of /| with decreasing T indicates 



a progressive decrease of the localization length [17 1 
However, in contrast to the standard behavior in crys- 
tal states, and similarly to relaxation in supercooled 
fluids, hopping events move the particles beyond the 
localization length and ergodicity for self-motions is 
restored, yielding final relaxation of Fs{q,t). 

The inset in Fig. [Ifb) shows that, in the same T- 
range of the latter data, the plateau height, fq, for 
F{q, t) stays constant, contrary to the case of /|. This 
behavior is observed for all q- values, except for those 
probing the lattice structure, for which a progressive 
increase of the plateau is observed with decreasing T 
(not shown). An increase of fq for q not probing the 
lattice structure is only observed at much lower T [see 
data for T = 0.070 in Fig. [ijb)]. Though for these 
low-T values the system cannot be equilibrated and 
aging effects are observed at long times (not shown), 
such effects do not affect the behavior of F{q,t) at 
intermediate times, where the plateau regime arises. 
Results in Fig. EJb) show that the localization tran- 
sition occurs at lower T for out-of-lattice collective 
correlations than for self-correlations, i.e., there is a 
dynamic decoupling between self- and collective re- 
laxation. These results resemble dynamic features of 
plastic crystals, where molecules are constrained to 



vibrate around lattice positions, but can perform full 
rotations leading to relaxation of out-of-lattice collec- 
tive correlations (in the present case cluster deforma- 
tion is an additional mechanism). However, contrary 
to the case of plastic crystals, activated hopping re- 
stores ergodicity of translational motions and leads to 
finite diffusivity and full decay of self-correlators. 

Results for fq (excluding g-values probing the lat- 
tice structure) above its localization transition (i.e., 
T-independent) are shown in Fig. [7] for the investi- 
gated densities. The approximate scaling behavior of 
fq is consistent with a common lattice constant for 
all the densities. Since ric oc p, the major effect of 
increasing Uc on F{q,t) is a rescaling of pq, which is 
canceled after normalization of F{q,t). The smaller 
width of fq as compared to /| indicates a much weaker 
localization for collective than for self-motions. 

In summary, we have investigated slow dynamics in 
the cluster crystal phase of a representative model of 
macromolecules in solution interacting through effec- 
tive ultrasoft bounded potentials. The obtained major 
features are a full change of the initial identity of the 
clusters through particle hopping between lattice sites, 
and dynamic decoupling between self- and collective 
out-of-lattice correlations. 
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